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Recent full-particle simulations of electron-positron reconnection have revealed that the Weibel 
instability plays an active role in controlling the dynamics of the current layer and maintaining 
fast reconnection. A four-beam model is developed to explore the development of the instability 
within a narrow current layer characteristic of reconnection. The problem is reduced to two coupled 
second-order differential equations, whose growing eigenmodes are obtained via both asymptotic 
approximations and finite difference methods. Full particle simulations confirm the linear theory 
and help probe the nonlinear development of the instability. The current layer broadening in the 
reconnection outflow jet is linked to the scattering of high- velocity streaming particles in the Weibel- 
generated, out-of-plane magnetic field. 



I. INTRODUCTION 

The temperature anisotropy driven Weibel instability 
[T is thought to play an important role in several as- 
trophysical systems. For instance, Weibel-mediated col- 
lisionless shocks in relativistic jets, pulsar winds, and 
gamma-ray bursts have been suggested (|2] [3] [4]) as a 
possible particle acceleration mechanism. The Weibel- 
generated magnetic field scatters particles, enabling them 
to bounce back and forth across the shock front, leading 
to acceleration via the first-order Fermi mechanism. [5 

Recently the role of the Weibel instability in electron- 
positron (pair) reconnection has begun to receive no- 
tice. Magnetic reconnection is a fundamental problem in 
plasma physics, and is ubiquitous in astrophysical phe- 
nomena involving magnetic fields, where it is the princi- 
pal mechanism for transforming magnetic energy into ki- 
netic energy and heat. Historically, the greatest difficulty 
in modeling reconnection has been to demonstrate that it 
is fast enough to match observations of energy release in, 
for instance, solar flares. By comparing multiple simula- 
tion models (e.g., two- fluid, hybrid, and full particle-in- 
cell (PIC)) the GEM challenge [6 demonstrated that in- 
clusion of the Hall term in the generalized Ohm's law was 
sufficient to produce fast reconnection. However, recent 
studies of electron-positron reconnection (in contrast to 
the usual electron-proton case) raised serious questions 
about the necessity of the Hall term for producing fast 
reconnection. In contrast with electron-proton plasmas, 
the mass symmetry in pair plasmas eliminates the Hall 
term. Yet, simulations suggest pair reconnection is still 
fast. Bessho & Bhattacharjee [7 attribute this fact to 
contributions from the off-diagonal components of the 
pressure tensor. Daughton & Karimabadi |8j later dis- 
cussed the role of island formation along the reconnec- 
tion layer. However Swisdak et al.[9 recently proposed 
that the Weibel instability, driven by an temperature 
anisotropy arising as inflowing plasma mixes with out- 
flow from the x-line, localizes the reconnection layer, and 
leads to fast pair reconnection. 

As an example of the importance of the Weibel in- 
stability in pair reconnection, we show, in Fig. [l] that 
suppressing the Weibel instability strongly influences the 



morphology of the current layer. In Fig. [TJa) we show a 
standard pair reconnection simulation where the Weibel 
instability causes the current layer to become turbulent 
and broaden, which opens the outflow exhaust as in 
Petschek's model [10]. In Fig.[ljb), we show the results 
of a simulation in which we suppress the instability by 
forcing the out-of-plane component of the magnetic field 
to zero. It is evident that without the turbulence pro- 
vided by the Weibel instability the narrow current layer 
extends to the system size. The longer current layer re- 
duces the reconnection rate in the simulation of Fig. [ijb) 
by one third. Since the Weibel instability plays such an 
important role in maintaining fast reconnection in pair 
plasma, a thorough understanding of its development in 
current layers is crucial. 

Although Swisdak et al. [9 proposed that the Weibel 
instability strongly influences reconnection, they only 
briefly considered the effects of the current layer envi- 
ronment on the instability's development. In pair re- 
connection the outflow layer is typically narrow (on the 
order of a few electron inertial lengths) and confined on 
both sides by regions of strong magnetic field. In this 
work we take a closer look at the role of the current layer 
and, in particular, how it slows (or prevents) the onset 
of the Weibel instability. We also try to understand how 
the unstable Weibel mode is able to open (broaden) the 
current layer. 

In section [TT] of this paper we introduce our analytic 
model and its assumptions. Section HI includes the 



derivation of the homogeneous dispersion relation and 
a comparison with kinetic theory. In Section IV we in- 



troduce the inhomogeneity arising from the reconnection 
geometry, and then numerically compute the eigenmodes 
and benchmark them with asymptotic theory in the lim- 
its of large and small current layer widths. In section \V\ 
we report on particle simulations that produce Weibel 
modes inside a Harris sheet. In section [Vl| the implica- 
tions for pair plasma reconnection are discussed and the 
downstream turbulent structure is compared with that 
of unstable Weibel modes. 
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II. THE GOVERNING EQUATIONS 

A cartoon of the basic Weibel instability is shown 
in Fig. [2j Consider a uniform unmagnetized plasma 
with beams counter-propagating in the x-direction (more 
generally, possessing an anisotropic temperature with 
Tx > Ty^Tz). If a sinusoidal magnetic field component 
Bz{y) arises from noise the positively charged particles 
with velocities Vx > ^ will converge towards the x-axis be- 
cause of the V X B Lorentz force while those with Vx <^ 
will diverge. The combination leads to a current density 
Jx{y) of the correct sign to amplify B^^ thus driving the 
mode unstable. (Negatively charged particles move in 
the opposite directions but have the same net effect.) 

To explore the structure of the Weibel instability in 
reconnect ion generated current layers, we use a fiuid de- 
scription. Since we are considering a pair plasma we in- 
clude four species, a G , e+, e— }, in our model: 
species of positrons and electrons with bulk velocities V 
and —V. In standard notation the governing equations 
are: 



d 



V-(n,Vc,) = 0; 



(1) 



m— Vc, = g«E+ ^Vc, X B ; (2) 

at c Tie 



V X B = ^ UaQa^a; 



(3) 



either electrons or positrons is 2n. If all perturbed vari- 
ables are proportional to e^^ we can linearize equations 
([l|-(|4| as 



jmVx.a + "rnVy^aVx^^ = qa^x] 



(7) 



(8) 



c n"^ n 

(9) 

imV^^a = QaE^ - '^Vy^aB^; (10) 



K = ^^qa{nv^,a + nv,,a); (n) 



iB, = cK, (12) 

where a tilde indicates a perturbed quantity. 

The effective x-direction temperature is Tx = tuV^ + 
Txx' Note that by assuming d/dx = we eliminate Txx 
from the equations. The effective temperature anisotropy 
is then determined only by the streaming velocity. 

Due to the mass and charge symmetry between elec- 
trons and positrons, we can collapse these eighteen equa- 
tions, Eq. ([7|)-Eq. (12), into two second-order differential 
equations (see Appendix A for details): 



dt 



B = -cV X E. 



(4) 



We assume the pressure tensor can be written in the di- 
agonal form 





" T 





o" 









T 





(5) 










T 





and that the temperature components do not vary in 
space or time. 

To particularize our coordinates we take the counter 
streaming velocities to be parallel to x. The perturbed 
magnetic field of the Weibel mode can then, without loss 
of generality, be taken to be parallel to z and the wavevec- 
tor parallel to y. All physical quantities are assumed uni- 
form in both the x and z directions, d/dz = d/dx = 0- 

Our initial state is characterized by 



V = F,(^)x; B = 5,(^)x; 

'^{y)'-) Tyy^Q^ Tzz,a T. 



(6) 



Pressure balance requires {4.nTyy + B^/Stt) = 0, where 
a prime stands for d/dy^ and the total number density for 



n 



4e 



VxuE'^ = T^x; 



Ex 



Aire 



(13) 



(14) 



with the following definitions: x = ^(K/,p+ — ^y,p- " 
Vy^e+^Vy.e-)] ^ IS the gyrofrcqucucy based on Bx] Cs = 
^JT jm is the sound speed; and d = ^y mc^ / Sirne^ = c/cjp 
is the skin depth with cOp the plasma frequency. 



III. THE WEIBEL INSTABILITY IN A 
HOMOGENEOUS PLASMA 

In a uniform magnetized plasma a dispersion relation 
can be found by combining Eqs. ( 13) and (14) and letting 
d/dy iky'. 



(15) 



Clearly while the streaming temperature (the first term 
of the right hand side) serves as the instability driver, 
thermal effects in the y-direction (the second term) and 
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the background magnetic field (the third term) stabilize 
the mode. 

In the even simpler case of a strongly anisotropic un- 
magnetized plasma, the kinetic growth rate of the Weibel 
instability is [9^, [11^: 



7 



2 + kM^ ' 



(16) 



of the resulting matrix. We use a grid size of 5/100 and 
domain size of 205, both of which are sufficient to en- 
sure covergence. Before describing the numerical results, 
however, we investigate the behavior of the equation an- 
alytically. 



kyde > 1, S/de < 1 



where Vth,x = \/Tx/m is the thermal velocity in the x 
direction. It is evident that this matches the first term of 
Eq. (15) with Vth,x ^x- The validity of using counter 



streaming cold plasma beams to analyze a single warm 
plasma was demonstrated by Davidson et al. [12 , who 
showed that this instability is not affected by the detailed 
shape of the plasma distribution function, but only by the 
effective temperature. 



In this limit the current layer thickness and the wave- 
length of the instability are much smaller than the iner- 
tial length d^- From general considerations we expect the 
mode to be harder to excite in such circumstances, unless 
is large. It is straightforward to combine Eqs. ( 13 ) and 



( 14 ) into a single second-order ordinary differential equa- 



tion by eliminating the second term in Eq. (14), which is 



small since kyd^ ^ 1. The result is 



IV. THE WEIBEL INSTABILITY INSIDE A 
HARRIS SHEET 

The results of Section IIIII were derived in the context 
of a homogeneous plasma. However for reconnection sim- 
ulations it is important to study how the development of 
the Weibel instability proceeds inside a narrow current 
layer. In this section we seek to understand whether, and 
to what degree, the inhomogeneities in the plasma den- 
sity and background magnetic field affect the instability. 

Our equilibrium is taken to be the usual Harris 
layer with a slight modification that incorporates an 
anisotropic plasma temperature 



= Bxo tanh(?//(5); Up^e = 2n = nosech^ {y /S); 



= {V^l - C's>ech\y/S) + Cl (18) 

where B^o^^o^ and Vxo are constants, S measures the 
width of current sheet and the subscripts p/e stand for 
positron/electron. We will refer to this setup as Profile 
A. The profiles are shown in Fig. [sj^a). 

For the sake of comparing with full particle pair re- 
connection simulations, all physics quantities are pre- 
sented in the same normalized units as those in Swisdak 
et al. [9 : the magnetic field to the asymptotic value 
of the reversed field B^o^ the density to no which is the 
value at the center of the current sheet minus a possi- 
ble uniform background density, velocities to the elec- 
tron Alfven speed Va^c^ lengths to the electron inertial 
length de = Y^mc^/47rnoe^, times to the inverse electron 
cyclotron frequency and temperatures to rrieV^ ^. 

We use T = 0.25meVA,e in the analysis of this section, 
which is the initial set up of the pair reconnection simu- 
lation of [9 , except for the latter's inclusion of a uniform 
background density of 0.2no. 

To find the modes of Eqs. (13) and ([T4|) we discretize 



the governing equations in the y direction (imposing zero 
derivative boundary conditions) and find the eigenvalues 



c -x' 

n 



2^ 

d2 



l^X- 



(19) 



For the parameter regime of interest, d/d^ ~ ^/^ce ^ 
0(1). Therefore for every term in Eq. (19) to be of the 
same order (except, perhaps, the termj^ the following 
scaling rules must apply: 



7 ^ Vx/d ^ kyCs > 1^ Cs/d. 



(20) 



After substituting for the functional form of the inhomo- 
geneity and changing variables to u = cosh.{y/d)x, we 
have 



2. // 



Q{y)u, 



(21) 



(17) where 



tanh^(y/<5) 



(22) 



and = Bxoe/mc. For V^q > Cgdl/2S'^, Q is negative 
near y = and positive at large y and Eq. (21 ) therefore 
has bounded solutions. 

We further simplify Q by Taylor expanding for small 
y/S^ and neglecting the tanh^ and (^q tanh^ terms 
in Eq. (22), which are small in the ordering given in 



Eq. (M, 



d^J^ 



y u : 



2T/?Q C'i 2 



(23) 



This equation can be solved in terms of Hermite polyno- 
mials with eigenvalues that give a maximum growth rate 
of 



7' 

/ max 



dj 



52 ■ 



(24) 
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Clearly the density shear term, C^^jb^ ^ which arises from 
the second term of Eq. (19), has the strongest stabiliz- 
ing effect. By letting ^rnax = 0, we obtain the critical 
temperature anisotropy required for the unstable mode 
in the small J/dg limit. 



3 + V5 (de_ 



where T^o {Tx)y=o rnV^Q. 

We plot 7 versus anisotropy in Fig. [4]^ a) for 5 = O.lde 
and in Fig. ^h) for S = 2de. In the S = O.lde case, the 
full equations (Eqs. (13),(14)) and the reduced equation 
(Eq. ([21])) result in the same curve because the wavenum- 
ber of the growing mode is large enough to validate the 
approximation. The analytical solution from Eq. (24) 
follows the correct trend and matches the numerical re- 
sults, particularly in the large anisotropy limit. We ob- 
serve that 7 is proportional to Vxo (i.e., Vxo ~ (T^o/T)^^'^ 
) as the scaling rule (Eq. ([20|) suggests. By compar- 
ing the growth rate of the unmagnetized homogeneous 
plasma from Eq. (15), it is also clear that the instability 



is severely suppressed by the inhomogeneity. In contrast, 
the 6 = 2de case results in a somewhat closer match be- 
cause of the increase in the inhomogeneity scale length. 



B. kyde < 1, S/de > 1 

Although in this limit the search for bounded modes is 
rather complicated, it is possible to gain some insight by 
expanding the homogeneous dispersion relation, Eq. ( 15 ) 
in the small limit: 



v3e. ( 1 



If we only keep terms of 0{ky ) we find that 



(26) 



(27) 



which is always positive and has no bounded modes. It 
is only by keeping the next term (i.e., 0{ky)) in the ex- 
pansion that bounded modes exist, a fact that will guide 
our treatment of the full system. 

We begin by neglecting the E'^ term in Eq. ( [l4| , solv- 
ing for Ex , and then substituting the result back into the 
equation, ultimately giving 



2 F 



Aire 



47re / cP 



(28) 



We then use this approximation in Eq. (13) in the small 



y/S limit and assume that we are close to marginal sta- 
bility. The resulting equation is: 



-i^x- 



(29) 



The parameter regime of interest to us is Vxo/VA,e ^ 
Cs/Va^c ^ ^o/^ce ~ 0(1). By requiring each term in 
this equation to be the same order we arrive at the or- 
dering 



1' - kliVio - C!); 



2 j2 



(^.'o-Cf)/^.'o«l- 



(30) 



(25) To proceed, we transform the equation to Fourier space 



dkl 



X- 



9 '^y 



(31) 



The quanti ty inside the square brackets has a maximum 
at kyo = ^ V^Q — C^/de Vxo . Since we are looking for the 
maximally growing mode, we Taylor expand this quantity 
in 5 = k.j — k.jci around s = 0: 



(52 ds'^ 



CD 



X- 
(32) 

The solution of this equation can again be written in the 
form of Hermite polynomials with a maximal eigenvalue 
of 



7' 

/ max 



2\2 



S 



(33) 



Without the second term, which arises from the back- 
ground magnetic field, the growth rate scales as j'^^x ^ 
{V^Q — Cg)ky (with ky given in Eq. (30)), which is essen- 
tially the same result as the unmagneuzed homogeneous 
relation, Eq. (15), in the kyde <^ 1 limit. From Eq. (33), 
we derive the marginal criterion in the large 6/de limit. 



Txo 



f 2V2nodl \ 
I CsS 



2/3 



(34) 



Finally we plot the threshold of marginal instability 
for different values of 5 in Fig. |5j The numerical solution 
of our model with Profile A is carefully benchmarked by 
these asymptotic theories in both the small and large S 
(or ky) limits. The (blue) dashed and (red) dot-dashed 
lines are discussed later. 



V. THE SMALL BOX PIC SIMULATION 

In order to confirm our linear theory and study the 
nonlinear development of the Weibel instability, we con- 
duct several simulations with the particle-in-cell (PIC) 
code p3d ^13^. The electromagnetic fields are defined 
on gridpoints and advanced in time with an explicit 
trapezoidal-leapfrog method using second-order spatial 
derivatives. The Lorentz equation of motion for each 
particle is evolved by a Boris algorithm where the veloc- 
ity V is accelerated by E for half a timestep, rotated by 
B, and accelerated by E for the final half timestep. To 
ensure that V • E = Airp a correction electric field is cal- 
culated by inverting Poisson's equation with a multigrid 
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algorithm. Although the code permits other choices, we 
work with fully periodic boundary conditions. 

We consider a system periodic in the x — y plane. 
The simulations presented here are two-dimensional, i.e., 
d/dz — 0. The initial equilibrium consists of a Harris 
current sheet superimposed on an ambient population of 
uniform density n^, 

Bx = Bx,h tanh(?//(5); Up^e = nhsech^{y/S) + n^, 

(35) 

where Bx^h^^h^^b are constants, subscript h stands for 
Harris sheet and 5 is the half width of the current sheet. 
We input the Harris plasma with an initial temperature 
Txx,h > Tyy = Tzz = T. The background plasma has an 
isotropic temperature T. Therefore, 

_ UbT ^nhsech^{y/6)Txx,h 



Ub + n/,sech {y/5) 



(36) 



where T and Txx,h are constants. This equilibrium is very 
similar to Profile A defined in Eqs. (17)-(n!8|) except for 



the inclusion of a background density. In the following 
we take the parameters rih = 0.8no,n5 = 0.2no,5^ ^ = 
O.SB'^Q^S = 2de in order to compare later to our full 
simulations of pair reconnection, and refer to the initial 
condition as Profile B (see Fig. [sj^b)). 

In order to prevent the simultaneous growth of the 
tearing mode [14 we consider a domain with Lx = ^de <C 
Ly = 50de. We let the Weibel mode grow from noise for 
different values of the temperature anisotropy, measure 
the growth rate, and compare the results with the theo- 
retical values from the full four-beam model in Fig.[6]^b). 
Kinetic effects only slightly reduce the growth rates of 
our four-beam results and the eigenfunction predicted by 
the model is in good agreement with the simulations (see 
Fig. [6]^ a)). The anisotropy threshold for the four-beam 
model with Profile B is shown as a (blue) dashed line in 
Fig. in 

The evolution of the temperature anisotropy and 
Weibel-generated from four small-box runs with dif- 
ferent initial conditions is shown in Fig. |6|c)-(d). The 
anisotropics (at ^ = 0) decrease in time to slightly above 
the marginal value ~ 2.3, while the amplitude of Bz si- 
multaneously rises. The increase of Bz scatters the hot 
streaming plasma, reducing the central anisotropy. Only 
a small part of the energy released transfers to Bz. The 
scattering increases Tyy^ and therefore the central pres- 
sure, causing the layer to expand. As a result, the ambi- 
ent Bx increases due to compression (not shown). 

In homogeneous plasmas we can predict the satura- 
tion level by analyzing the particle motion in the y- 
direction, dVy/dt = —{e/mc)VxBz. Roughly speaking, 
the Weibel instability saturates when the magnetic field 
grows to a value such that the particles become magnet- 
ically trapped and can no longer amplify the field. Trap- 
ping occurs when the particle excursion along ^, Ay, is 
comparable to the wavelength during the mode growth 



time. 



Equivalent ly. 



Ay 



eVxBz 



( 6 

^raax ^ ( ^xBzky raax 

\mc 



1/2 



(37) 



(38) 



The saturation occurs when the magnetic bounce fre- 
quency, cjb, is comparable to the fastest linear growth 
rate. This is essentially the empirical result from David- 
son et al.[T2]. In the strong anisotropy limit, ^max ^ 
(Eq. (15)), we obtain the simple saturation cri- 
1, where pz is the gyro-radius in the 



^xky,max 

terion ky 



^y^raaxPz 

Bz field. In our inhomogeneous plasmas the same prin- 
ciple guides the saturation of the Weibel mode, although 
the predicted saturation value of Bz is smaller due to the 
inhomogeneity. 



VI. IMPLICATIONS FOR PAIR 
RECONNECTION 

A. ky structure 

The basic features of a pair plasma reconnection sim- 
ulation are shown in Fig. [71 The Weibel instability 
manifests itself as a chess-board-like structure in the 
downstream out-of-plane magnetic field (see Fig. [7|(b)). 
The structure travels downstream with the outfiow speed 
from the x-line, which implies that the Weibel is a purely 
growing mode in the frame of the outfiowing plasma. 
Swisdak et al. [9] proposed that the anisotropy (with typi- 
cal magnitude Tx/Ty ~ 2—4) driving the instability arises 
from cold infiowing plasma mixing with outflow from the 
X-line. Here we apply our analytical results to these ob- 
servations. After fitting the inhomogeneity (F ig. jT^d)) 
seen in the reconnection simulations by Eqs. ([35^(36) 
with parameters Uh = 0.25no, rib = 0.2no, B^ ^ = O.SBJq 
(denoted Profile C) we can find the necessary temper- 
ature anisotropy for Weibel to be unstable within the 
current layer. The (red) dot-dashed curve in Fig. [s] is 
the marginal criterion for the Weibel instability, and it 
gives us the minimum temperature anisotropy Tx/T ~ 
3.3 for a typical reconnection layer of width Me (i.e., 
6 = 2de). This predicted minimum is within the observed 
anisotropy values (2 ~ 4) and thus demonstrates that 
the reconnection current layer can support the Weibel 
mode. The even parity of Bz produced by the Weibel 
mode (Fig. |6|a)) is also seen in the reconnection sim- 
ulation. Farther downstream, the nonlinear saturation 
of the Weibel instability stops the system from moving 
to even higher temperature anisotropics and keeps the 
system near marginal stability. 

The magnitude of Bz saturates at ~ 0.1 — {)ABxq in 
Fig. j6|d), which is comparable to the value shown in 
Fig. ufb). Therefore by measuring the time, t^, for Bz 
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to saturate, we can estimate the half-length of the re- 
connection nozzle as L ^ VA,ets since Bz develops from 
the instability while the plasma is convected out from 
the x-point at roughly the electron Alfven velocity. As a 
result, anisotropy values of 2.5 ~ 4 with saturation times 
from Fig. |6jd) give a predicted half-nozzle length in the 
range 25 ~ SOd^^ which compares favorably with the ob- 
served reconnection layer half-length of GOdg- Although 
these small runs use higher density plasmas than the re- 
connection runs, we expect the saturation behavior to 
be similar. Furthermore, the (red) dot-dashed curve in 
Fig.|6](b), which represents the four-beam model solution 
with Profile C, indicates that an anisotropy of 4.0 pro- 
duces a growth rate of ~ 0.2Qce^ which again leads to a 
nozzle half-length of ~ 50de (i.e., we expect its evolution 
to be similar to the solid curve in Fig. |6](d)). 

In contrast to our small-box runs where reconnection 
was suppressed, we expect the background to be nois- 
ier in a simulation that allows both reconnection and 
the Weibel instability to develop. However that will not 
strongly affect our results. Since the out-of-plane field Bz 
grows exponentially from the initial noise we expect the 
estimated nozzle length to scale logarithmically with the 
initial noise level. Hence the predicted length via this es- 
timate is insensitive to the noise in Bz and we can expect 
our small-box runs to still give a reliable comparison to 
reconnection simulations. 



B. kx structure 

We now focus on the origin and magnitude of the finite 
kx in the region downstream of the x-line. Loosening our 
assumptions by letting d/dx 7^ in the linearized four- 
beam model from Eq. Q-Q, we can numerically solve 
for the two-dimensional dispersion relation in a homoge- 
neous plasma. One new feature is the introduction of a 
two-stream instability (i.e., ky = mode) that has previ- 
ously been noted by Zenitani & Hesse [15 (see Appendix 
B for the dispersion relation). As shown in Fig. [sj the 
two-stream instability has a higher growth rate than the 
Weibel instability. Hence it always grows in front of the 
Weibel structures with wavelength ~ 27r/3 ~ 2de^ where 
the factor of 3 comes from the wavenumber for the max- 
imum growth mode measured in Fig. Sj^c). We do ob- 
serve a double peaked distribution of the x-direction ve- 
locity and the signature of the two-stream instability 
(Fig. [t] (e), (c)) in the appropriate region with a wave- 
length comparable to the predicted value. The nonlinear 
development of the two-stream instability tends to merge 
the counter-streaming distributions, resulting in a single- 
humped distribution farther downstream with > Ty. 
The result is a transition from coexisting Weibel and two- 
stream instabilities to a pure Weibel instability, as can 
be seen in Fig. ^h). 

The two-stream maximum growth rate predicted by 
the four-beam model is three times larger than that of 
the Weibel instability (divide the maximum growth rate 



of the solid curve in Fig.|8jc) by that of Fig.jsjb)), which 
seems inconsistent with their relatively close development 
in the reconnection simulations. However this disagree- 
ment can be reduced by the introduction of kinetic effects 
and a finite T^x^ both of which are not included in our 
four-beam model. In particular, a finite Txx suppresses 
the growth of the two-stream instability. Kinetic theo- 
ries predict lower growth rates for both the two-stream 
and Weibel instabilities, with the ratio of their fastest 
growth rates decreasing to about 1.6 (divide the maxi- 
mum growth rate of the dashed curve in Fig.jsjc) by that 
of Fig.jsj^b)). This lower ratio helps to explain the nearly 
coincident signatures of both instabilities in pair recon- 
nection. In general, we also expect that a full kinetic 
treatment of Harris sheet inhomogeneous plasmas will 
produce slightly lower Weibel growth rates than those of 
our four-beam model, a feature that has already been 
observed in Fig. [6]^b). 

The two-stream instability can not explain the longer 
x-direction variance of the chess-board-like structure far- 
ther downstream. In that region, Bx decreases from the 
asymptotic value of the reversed field Bxo) at the 
nozzle edge to zero at the nozzle symmetry line while the 
temperature anisotropy remains large. Although Swis- 
dak et al. [9^ use a parity argument to argue against 
the possibility that the firehose instability plays a role, it 
is perhaps plausible that it could couple to the stronger 
Weibel mode and provide the finite kx. However this pos- 
sibility is again ruled out by directly solving for the ho- 
mogeneous dispersion relation of the firehose instability 
via a one fluid double-adiabatic model with finite Larmor 
radius corrections [16 . The firehose instability in a ho- 
mogenous plasma with Up^e = 0.45no, Tx = 4.0T has the 
strongest growth rate ~ O.OSb^ce at Bx ~ 0.445^^0 and 
kxde ~ 0.07, which is an order of magnitude smaller than 
the growth rates of the Weibel and two-stream instabil- 
ities. Moreover, the predicted wavelength is too large to 
explain the observed kx- 

A possible mechanism for the observed kx is proposed 
here. First imagine that the instability is confined in the 
y direction between reflecting walls (see Fig. [9|. The re- 
flecting walls mimic the strong Bx field at the bound- 
aries of the current layer. In this geometry there is 
an intrinsic scale length A^^ associated with the trajec- 
tory of a particle away from the current layer, its reflec- 
tion from the wall, and its motio n back towards current 
layer. For a uniform Bz^ Ax = p\ — [pz — A)^, where 
pz = Vx/^c,z is the gyro-radius based on Bz and A is the 
half distance between the walls. The particle trajectory 
is sketched in the dashed box of Fig. [9ja). We suggest 
that this intrinsic scale length controls the x-dependence 
of Bz seen in the simulations. We define the scale A^^ 
shown as a dashed box in Fig. [9]^a) as a Weibel-unit and 
consider the interaction between two Weibel- units. When 
all of the converging streaming particles (dashed curves 
with rightward arrows) leave the center of the righthand 
Weibel-unit in Fig. [9]^ a), this unit will need a replenish- 
ment of converging streaming particles (solid curves with 
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rightward arrows) from its leftward neighbor to maintain 
its central current. Hence this self-consistent arrange- 
ment can arise without x-variation. 

However, suppose that due to the initial random noise 
the lefthand Weibel-unit acquires the opposite polarity 
magnetic field, as seen in Fig. |9jb). Then the intrinsic 
will arise such that the source for replenishing the 
converging streaming particles of the righthand Weibel- 
unit is the diverging streaming particles (solid curves 
with rightward arrows) of a leftward neighbor after they 
bounce (perhaps multiple times) against the walls. The 
lefthand Weibel-unit can be reinforced from its right- 
ward neighbor in the same manner. In this configuration 
we can estimate the distance between two neighboring 
Weibel-units to be NA^^ where TV is the number of times 
a particle refiects from the walls. 

This qualitative explanation is supported by PIC sim- 
ulations, as can be seen in Fig. 



10 We confine an 



plasma with an anisotropic temperature within a mag- 
netic trough with thickness 5de. Specifically, Bx = 0^ 
Tx 4.0T, and n 1.2no between y -2.5 ^ 2.54 
and Bx = l.OBxo^ Tx = T, and n = 0.2no outside 
this region. The simulation is performed in a domain 
of size 200de x 25de with periodic boundary conditions 
in both the x and y directions. Here the Bx trough 
serves as the refiecting walls and the trough thickness 
of 5de is comparable to the reconnection nozzle thick- 
ness. Two distant Weibel-units are not able to commu- 
nicate with each other, and therefore it is not surpris- 
ing to see that the magnetic fields at x = —bOde and 
bOde in Fig. [lOja) have opposite signs. As time evolves, 
the Weibel-generated field gets stronger while the long 
Weibel-unit {x = — SOdg ~ 30de at Fig. lOjb)) breaks up 
into smaller Weibel units of opposite polarity at (c). Five 
test positrons with initial velocity y^^TlJm = I.OVa^c 
(i.e., Tx/T = 4) are randomly placed near (lOdg, Ode). 



Their trajectories are shown as white curves in Fig. 10 'c) 
and are blown up in Fig. 10 'd). They are qualita- 
tively similar to those described in Fig. |9jb). The con- 
verging motions of those trajectories at the trough cen- 
ter between x = 10 ^ SOd^ and x = 45 ~ 65de of 



Fig. 10 ^d) contribute to the current for the out-of-plane 
magnetic field. These particles refiect from Bx in be- 
tween X = 30 ~ 45(ie, where they deposit their momenta. 
At late time, tQce = 37.5, the entire channel relaxes to 
a chain of Weibel units of alternating polarity, whose 
length scale is approximately determined by A^A^^ with 

In order to apply this idea to the downstream regions 
in reconnection simulations, we approximate the charac- 
teristic streaming velocity as 1.0yA,e and the averaged 
out-of-plane magnetic field as 0.2^3^0- This implies a 
gyro-radius of bd^. The layer thickness A is 2de along 
the reconnectio n nozzle, and thus A^; has a scale of ap- 
proximately 2y^5^ ~ (5 — 2)^ = 8(ie, about half of the 
size of the observed structure in Fig. [tJ^c). This mecha- 
nism roughly explains the scale size of the observed vari- 
ation. Furthermore, the gyromotion of the reconnection 



outfiow helps explain how the reconnection nozzle broad- 
ens downstream. The out-of-plane magnetic field bends 
the outfiow momentum from the x to the ^/-direction. As 
a result this fiow pushes the background Harris magnetic 
field away from the symmetry line, broadening the cur- 
rent layer. 



VII. SUMMARY AND DISCUSSION 

We have developed a four-beam model to study the ef- 
fects on the Weibel instability of a spatial inhomogeneity 
arising from a current layer. We have shown that the 
Weibel instability is able to grow within narrow Harris 
sheets, and its growth rate (saturation time), saturation 
magnitude, and mode structure fit those values observed 
in pair reconnection. 

This further suggests that the Weibel instability might 
control the current layer dynamics in pair reconnection, 
where the Hall term is absent. Other candidate insta- 
bilities, such as the two-stream and firehose instabilities 
have rather minor effects, particularly since it is not clear 
if the firehose instability even appears in these systems. 
The high- velocity outfiow scatters into the transverse di- 
rection due to the Lorentz force arising from the Weibel- 
generated out-of-plane magnetic field. We argue that the 
associated increased pressure Pyy is responsible for open- 
ing the pair reconnection nozzle and shortening the cur- 
rent layer. As a consequence, the shorter current layer 
generates a higher reconnection rate. 

Even though a similar temperature anisotropy also 
could arise in electron-proton reconnection, signatures of 
the Weibel instability are not seen there. The reason is 
that the even shorter current layer (~ lOdg), controlled 
by whistler waves [17 , leaves insufficient space for the un- 
stable Weibel mode to grow. It is an open question as to 
how the Weibel instability controlled reconnection trans- 
forms into whistler mediated reconnection as the electron 
to ion mass ratio changes. 

The development of the Weibel instability in the ini- 
tial state of relativistic {ksT ^ mc^) pair reconnection 
is described in Zenitani & Hesse [15]. There the Weibel 
instability is shown to grow in front of a tangential dis- 
continuity formed by the mixing of outffowing and am- 
bient plasma, but it is not clear if this turbulence plays 
a role in the steady-state development of the outffow ex- 
haust and the associated current layer. However, by ar- 
tificially suppressing the Weibel-generated out-of-plane 
magnetic field, they do demonstrate the ability of the 
Weibel mode to broaden the current layer, similar to the 
behavior of the present non-relativistic case. The growth 
rate and wave vector of the relativistic Weibel instability 

1 /2 

are smaller by a factor of 7^ , where 7l is the Lorentz 
factor [18 . Consequently, we expect the instability to 
grow more slowly during relativistic pair reconnection, 
not only because of the intrinsically lower growth rate 
but also because of the relatively larger suppressing ef- 
fect of the Harris reversed field on the enlarged mode 
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Profile A 



(b) Profile B 



FIG. 1: (Color online) Out-of-plane electron velocity for two 
pair reconnection simulations. In (a) the current layer length 
is about ^ 120de (i.e., |x| < 60de) [2- In (b) the current layer 
length of a simulation with Bz = scales as the system size, 
^ 700de (i.e., |x| < 350de). The structure in the large island, 
which may be due to a two-stream instability, does not affect 
the behavior near the x-line. 
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FIG. 2: Schematic of the Weibel instability. 



FIG. 3: (Color online) Current sheet profiles. (Blue) solid 
curves are B^, (green) dashed curves are the positron/electron 
density np,e, and (red) dot-dashed curves are x-direction tem- 
perature Tx. Profile A in (a) with no = l.Ono, Bxo = l.OBxo, 
and 6 = 2de. Profile B in (b) with rih = 0.8no,n6 = 
0.2no, Bl^h = O.SBlo and S = 2de. 



(a) 5 = 0.1cl 




14 1 



(b) 8 = 2d 



■ Inhomogeneous 
Homogeneous 




structure (if we assume a similar nozzle thickness in both 
the relativistic and non-relativistic regimes). Overall, it 
remains an open question as to whether the Weibel in- 
stability plays an important role in controlling relativistic 
electron-positron reconnection. 



FIG. 4: Maximum growth rate versus temperature anisotropy 
for Profile A for (a) S = O.lde and (b) S = 2de. Thick 
solid curves: numerical solutions of the full four-beam model, 
Eqs. (13)- (14). Thin solid curve (only in (a)): numerical so- 
lutions of the reduced equation, Eq. (19). Diamonds (only 
in (a)): analytical solutions in the large-k limit, Eq. (24). 



Dashed curves: analytical solutions for an unmagnetized ho- 
mogenous plasma, Eq. (15). 
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VIII. APPENDIX A: THE DERIVATION OF 



THE GOVERNING EQUATIONS (13) AND (14) 



To derive Eqs. (13) and (14), we define the variables, 

(39) 



X^Vy 



(40) 
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Marginal criterion 



— Profile A 

- - Profile B 
- Profile C 




log(6/d^) 



FIG. 5: (Color online) The marginal threshold of the Weibel 
instability as a function of S. Profile A: the solid curve from 



the numerical solutions of the full four-beam model (Eqs. ( 13 )- 
( 14 )), the diamonds from the analytical solution in the small-^ 
limit (Eq. (25)), an d st ars from the analytical solution in the 
large-^ limit(Eq. (34)). Profile B: the (blue) dashed curve 
from numerical solutions of the full four-beam model. Profile 
C: the (red) dot-dashed curve from numerical solutions of the 
full four- beam model. 




FIG. 6: (Color online) Results from PIC simulations with 
initial Profile B. In (a) the eigenfunction of Bz with Txo/T — 
3.0: in solid the small-box PIC simulation and in (blue) dash 
the full four-beam model. In (b) the maximum growth rate 
versus anisotropy: in diamond the small-box PIC simulations, 
in (blue) dash the full four-beam model and in (red) dot-dash 
the full four-beam model with Profile C. In (c) the temporal 
evolution of different initial temperature anisotropics in the 
small-box PIC simulation. In (d) the magnitude of Bz (at 
2/ = 0) for the temperature anisotropics plotted in (c). 



V, 



z,e+ 



(41) 



T] = np+ - rip- 



■ + Tie 



(42) 



then combine Eqs. ([7|)-([T2| to yield 




' '20 
„ 10 




(43) 
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-77) = n'x + ^x', 


-20 





.6 (d) 1.0 



7mC 4.eEx - mV^x, 



p ~ 4p ~ n' T 

C C ~- 



n 



c 



~ / 47re 47re ~ 

Bz VxT] H nC, 

c c 



jB^ cE^ . 



(44) 
(45) 
(46) 

(47) 
(48) 



Eqs. (|43|)-([48| represent a set of 6 equations for the 6 vari- 
ables, {fjX^x^^^ Bz^ Ex}. Note that {n^Vx^ B^^T} are 
unperturbed quantities specifying the initial conditions. 
Now use Eqs. (|43|), (|44|), and (psl) to rewrite Eq. (|45|) in 



terms of x and Ex 



7^X ^ 



m 



(49) 





(c) 
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FIG. 7: (Color online) A PIC simulation of pair reconnection. 
In (a) the x-direction positron temperature. The X-point is 
on the right edge of the plot. In (b) the Bz signatures of 
both the Weibel (chess-board-like structure) and two-stream 
(finer structure upstream) instabilities. In (c) the Ex sig- 
nature of the two-stream instability in the downstream re- 
gion. In (d) the inhomogeneity plotted along the white line 
in (a). The line styles (colors) are the same as in Fig. p] In 
(e) the double-humped velocity distribution function (blue) 
in the two-stream and Weibel unstable region (x/de,y /de G 
(—21 : 18, —4 : 4)) becomes single-humped (black) farther 
downstream {x/de,y/de G (-44 : -30, -4 : 4)). 




(a) 2-D Dispersion relation 




(b) Weibel instability (k =0) 



(c) Two-stream instability (k =0) 





FIG. 



Homogeneous plasma dispersion relation with c 



0, = 4T^, 



-- 4.0meVA,e- In (a) 
the growth rate as a function of kx and ky . In (b) the homoge- 
neous dispersion relation for the Weibel instability. The solid 
curve corresponds to the four-beam model and the dashes to 
kinetic theory. In (c) the homogeneous dispersion relation 
of the two-stream instability. The solid curve corresponds to 
a the four-beam model with Txx = 0. The dashed curve is 
kinetic theory with Txx = Tyy = 0.25meVA e- 
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FIG. 10: (Color online) Evolution of the Weibel instability 
with initial anisotropy 4.0 inside a magnetic field trough. 
In (a)-(c), Bz at tQce = 12.5,25,37.5. White curves mark 
the trajectories of five test positrons with initial velocity 
yjTxjvn — 1.0VA,e and initial position (lOde, Ode)- In (d) 
enlargement of the particle trajectories at time 37.5/r2ce. 



Weibel-unit 




FIG. 9: Schematic explanation of the interaction between 
Weibel-unit s. White regions have > 0, hatched Bz < 0. A 
Weibel-unit is represented by particle trajectories similar to 
the dashed curves of Fig. (2] In (a) the converging particles of 
the right Weibel-unit are replenished by the converging parti- 
cles from the left unit. In (b) converging particles of the right 
Weibel-unit can be replenished by diverging particles from the 
left unit of opposite polarity. 



Use Eqs. (43), (46), and (48) to rewrite Eq. (47) in 



terms of x E^^ 



eJ' = —VAnxY + ^E. - ^V^nx. (50) 



After some minor algebraic manipulations, Eqs. (49) 
and ([sol ) can be written in the form shown in Eq. (13) 
and (fill. 



IX. APPENDIX B: THE DISPERSION 
RELATIONS OF THE TWO-STREAM 
INSTABILITY AND THE WEIBEL INSTABILITY 



If we consider the ky = O^kx ^ limit of our four- 
beam model, we arrive at the usual two-stream instability 
dispersion relation. 



(51) 



where co = cOr -\- ij . 

The kinetic version of this dispersion relation with fi- 
nite Txx is 

fc^C^ + c.2[2 + 6^(6) + 6^(6)] = 0, (52) 

where Z{^) = (l/^/n) J^^exp{-x'^)/{x - i)dx, 
a = {ujl\k^\ - V^)/^2Tj^/m and 6 = {'^/\kx\ + 
V;)/V2T^Jm. 

For reference, the dispersion relation for the Weibel 
instability {kx = 0,ky ^ 0) in kinetic theory is 



ky -u,^ + 2u>i[i ^ 



2loI^^Z{£), (53) 



where ^ = ij / {\ky\^2Ty / m) . Note that this reduces to 
Eq. (16) in the strong anisotropy limit. 
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If both kx and ky ^ 0, both instabihties are present. We treat this hmit numericaUy because of the complexity. 
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